Library loading¶

In [1]:
import scanpy as sc, anndata as ad, numpy as np, pandas as pd
from scipy import sparse
from anndata import AnnData
import warnings
import socket
import plotly.express as px
from matplotlib import pylab
import sys
import yaml
import os
import matplotlib.pyplot as plt
import scvelo as scv
from pybiomart import Server
import rpy2
from pandas.api.types import CategoricalDtype
import rpy2.robjects as ro
import seaborn as sns
from matplotlib_venn import venn2, venn2_circles, venn2_unweighted
from matplotlib_venn import venn3, venn3_circles
from matplotlib import pyplot as plt
import plotly.graph_objects as go
import ipyparams


warnings.filterwarnings('ignore')
/group/testa/Users/davide.castaldi/Polaroid_publication_revisions/.local/lib/python3.8/site-packages/setuptools_scm/_integration/setuptools.py:30: RuntimeWarning: 
ERROR: setuptools==56.1.0 is used in combination with setuptools_scm>=8.x

Your build configuration is incomplete and previously worked by accident!
setuptools_scm requires setuptools>=61

Suggested workaround if applicable:
 - migrating from the deprecated setup_requires mechanism to pep517/518
   and using a pyproject.toml to declare build dependencies
   which are reliably pre-installed before running the build tools

  warnings.warn(
In [2]:
outdir="./outdir"
FinaLeaf = "/Cycling"
output_basename = "09C_Cycling_pBulk_bySegment.DEA"
mappingDict = {}
categoriesOrder = ["piece1","piece2","piece3","distal","medial","proximal"]

groupingCovariate = "group"
PseudooReplicates_per_group = 10



totalPath = outdir+FinaLeaf+"/5C_Cycling_pBulk.bySegment."+str(PseudooReplicates_per_group)+"PRs.tsv"
adataPath = outdir+FinaLeaf+"/5C_Cycling_pBulk.bySegment."+str(PseudooReplicates_per_group)+"PRs.h5ad"

if not os.path.exists('./figures'):
    os.makedirs('./figures')

if not os.path.exists('./tables'):
    os.makedirs('./tables')
In [3]:
%load_ext rpy2.ipython
In [4]:
%%R
library(edgeR)
library(org.Hs.eg.db)
library(AnnotationDbi)
library(stats)
library(topGO)
R[write to console]: Loading required package: limma

R[write to console]: Loading required package: AnnotationDbi

R[write to console]: Loading required package: stats4

R[write to console]: Loading required package: BiocGenerics

R[write to console]: Loading required package: parallel

R[write to console]: 
Attaching package: ‘BiocGenerics’


R[write to console]: The following objects are masked from ‘package:parallel’:

    clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
    clusterExport, clusterMap, parApply, parCapply, parLapply,
    parLapplyLB, parRapply, parSapply, parSapplyLB


R[write to console]: The following object is masked from ‘package:limma’:

    plotMA


R[write to console]: The following objects are masked from ‘package:stats’:

    IQR, mad, sd, var, xtabs


R[write to console]: The following objects are masked from ‘package:base’:

    anyDuplicated, append, as.data.frame, basename, cbind, colnames,
    dirname, do.call, duplicated, eval, evalq, Filter, Find, get, grep,
    grepl, intersect, is.unsorted, lapply, Map, mapply, match, mget,
    order, paste, pmax, pmax.int, pmin, pmin.int, Position, rank,
    rbind, Reduce, rownames, sapply, setdiff, sort, table, tapply,
    union, unique, unsplit, which.max, which.min


R[write to console]: Loading required package: Biobase

R[write to console]: Welcome to Bioconductor

    Vignettes contain introductory material; view with
    'browseVignettes()'. To cite Bioconductor, see
    'citation("Biobase")', and for packages 'citation("pkgname")'.


R[write to console]: Loading required package: IRanges

R[write to console]: Loading required package: S4Vectors

R[write to console]: 
Attaching package: ‘S4Vectors’


R[write to console]: The following object is masked from ‘package:base’:

    expand.grid


R[write to console]: 

R[write to console]: Loading required package: graph

R[write to console]: Loading required package: GO.db

R[write to console]: 

R[write to console]: Loading required package: SparseM

R[write to console]: 
Attaching package: ‘SparseM’


R[write to console]: The following object is masked from ‘package:base’:

    backsolve


R[write to console]: 
groupGOTerms: 	GOBPTerm, GOMFTerm, GOCCTerm environments built.

R[write to console]: 
Attaching package: ‘topGO’


R[write to console]: The following object is masked from ‘package:IRanges’:

    members


In [5]:
%matplotlib inline
In [6]:
sc.settings.verbosity = 3         # verbosity: errors (0), warnings (1), info (2), hints (3)
sc.logging.print_header()
sc.settings.set_figure_params(dpi=50, facecolor='white')
pylab.rcParams['figure.figsize'] = (10, 10)
scanpy==1.8.0 anndata==0.8.0 umap==0.4.6 numpy==1.22.2 scipy==1.10.1 pandas==1.2.3 scikit-learn==0.24.2 statsmodels==0.13.1 python-igraph==0.9.8 louvain==0.7.1 pynndescent==0.5.5

Configure paths¶

In [7]:
hostRoot = "-".join(socket.gethostname().split('-')[0:2])



#indir=paths["paths"]["indir"][hostRoot]
#projectBaseDir=paths["paths"]["projectBaseDir"][hostRoot]

Load ExcM neurons pseudobulk and set categories order¶

In [8]:
adata = sc.read_h5ad(adataPath)
adata.obs[groupingCovariate] = adata.obs[groupingCovariate].astype(CategoricalDtype(categories=categoriesOrder, ordered=True))
In [9]:
total = pd.read_csv(totalPath, sep="\t", index_col=0)

Reorder counts cols¶

In [10]:
#adata = adata[adata.obs["group"].isin(RelevantAreas)]
total = total[adata.obs_names]

DEA¶

RPY INITIALIZATION¶

In [11]:
edgeR_topTags = ro.r['topTags']
edgeR_glmLRT = ro.r['glmLRT']
edgeR_glmQLFTest = ro.r['glmQLFTest']
as_data_frame = ro.r['as.data.frame']
rprint = rpy2.robjects.globalenv.find("print")
In [12]:
DEGs = {}
TestCov = groupingCovariate
In [13]:
adata.obs.group.value_counts()
Out[13]:
piece1      10
piece2      10
piece3      10
distal      10
medial      10
proximal    10
Name: group, dtype: int64

MM¶

In [14]:
obs = adata.obs[[TestCov,"pseudoreplicate"]].copy()
obs["pseudoreplicate"] = obs.index.tolist()
obs["TestCov"] = obs[TestCov].astype(str)

#obs["TestCov"] = obs[TestCov]

Counts¶

In [15]:
totalRelevant = total[obs.index]
totalRelevant.head()
Out[15]:
distal_0 distal_1 distal_2 distal_3 distal_4 distal_5 distal_6 distal_7 distal_8 distal_9 ... proximal_0 proximal_1 proximal_2 proximal_3 proximal_4 proximal_5 proximal_6 proximal_7 proximal_8 proximal_9
MIR1302-2HG 0.000000 0.0 0.000000 0.0 0.000000 0.0 0.0 0.0 0.000000 0.0 ... 0.0 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
FAM138A 0.000000 0.0 0.000000 0.0 0.000000 0.0 0.0 0.0 0.000000 0.0 ... 0.0 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
OR4F5 0.000000 0.0 0.000000 0.0 0.000000 0.0 0.0 0.0 0.000000 0.0 ... 0.0 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
AL627309.1 1.909059 0.0 1.303118 0.0 1.871744 0.0 0.0 0.0 0.475036 0.0 ... 0.0 0.724288 3.736357 1.517416 0.611369 0.897955 0.747047 0.500132 2.488713 1.779022
AL627309.3 0.000000 0.0 0.000000 0.0 0.000000 0.0 0.0 0.0 0.000000 0.0 ... 0.0 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000

5 rows × 60 columns

In [16]:
minCounts = totalRelevant[totalRelevant > 0].min(axis = 1).min()

#Universe 1 count in at least 5% of samples per group
groupThreshold = round(pd.get_dummies(obs[TestCov]).T.sum(axis = 1) *0.05)
bolMatrix = (np.matrix(np.dot(pd.get_dummies(obs[TestCov]).T, (totalRelevant > 0).astype(int).T)) - np.matrix(groupThreshold).T > 0)
bolVect = (bolMatrix.sum(axis = 0) >= 1).A1
totalRelevant = totalRelevant.loc[bolVect]
len(totalRelevant)
Out[16]:
25569
In [17]:
universe = totalRelevant.index.tolist()
In [18]:
levelsToMap = adata.obs[groupingCovariate].cat.categories.tolist()
levelsToMap
Out[18]:
['piece1', 'piece2', 'piece3', 'distal', 'medial', 'proximal']
In [19]:
%%R -i obs -i totalRelevant -i levelsToMap -o dds 

mmVector<- factor(obs$TestCov, levels = c(levelsToMap))

mm <- model.matrix(~mmVector )
row.names(mm) <- colnames(totalRelevant)

dds <- DGEList(totalRelevant, group = obs$TestCov, genes = rownames(totalRelevant))

#Ennotate with entrez genes
anno <- select(org.Hs.eg.db, keys=rownames(dds$counts),columns=c("ENTREZID","SYMBOL"),keytype="SYMBOL")
colnames(anno) <- c("genes","entrez")
anno <- anno[!duplicated(anno$genes, fromLast=T), ]
dds$genes$entrez <- merge(x = data.frame(dds$genes), y = anno, by = "genes", all.x = TRUE)$entrez

dim(dds)
R[write to console]: 'select()' returned 1:many mapping between keys and columns

[1] 25569    60
In [20]:
%%R

mm
           (Intercept) mmVectorpiece2 mmVectorpiece3 mmVectordistal
distal_0             1              0              0              1
distal_1             1              0              0              1
distal_2             1              0              0              1
distal_3             1              0              0              1
distal_4             1              0              0              1
distal_5             1              0              0              1
distal_6             1              0              0              1
distal_7             1              0              0              1
distal_8             1              0              0              1
distal_9             1              0              0              1
medial_0             1              0              0              0
medial_1             1              0              0              0
medial_2             1              0              0              0
medial_3             1              0              0              0
medial_4             1              0              0              0
medial_5             1              0              0              0
medial_6             1              0              0              0
medial_7             1              0              0              0
medial_8             1              0              0              0
medial_9             1              0              0              0
piece1_0             1              0              0              0
piece1_1             1              0              0              0
piece1_2             1              0              0              0
piece1_3             1              0              0              0
piece1_4             1              0              0              0
piece1_5             1              0              0              0
piece1_6             1              0              0              0
piece1_7             1              0              0              0
piece1_8             1              0              0              0
piece1_9             1              0              0              0
piece2_0             1              1              0              0
piece2_1             1              1              0              0
piece2_2             1              1              0              0
piece2_3             1              1              0              0
piece2_4             1              1              0              0
piece2_5             1              1              0              0
piece2_6             1              1              0              0
piece2_7             1              1              0              0
piece2_8             1              1              0              0
piece2_9             1              1              0              0
piece3_0             1              0              1              0
piece3_1             1              0              1              0
piece3_2             1              0              1              0
piece3_3             1              0              1              0
piece3_4             1              0              1              0
piece3_5             1              0              1              0
piece3_6             1              0              1              0
piece3_7             1              0              1              0
piece3_8             1              0              1              0
piece3_9             1              0              1              0
proximal_0           1              0              0              0
proximal_1           1              0              0              0
proximal_2           1              0              0              0
proximal_3           1              0              0              0
proximal_4           1              0              0              0
proximal_5           1              0              0              0
proximal_6           1              0              0              0
proximal_7           1              0              0              0
proximal_8           1              0              0              0
proximal_9           1              0              0              0
           mmVectormedial mmVectorproximal
distal_0                0                0
distal_1                0                0
distal_2                0                0
distal_3                0                0
distal_4                0                0
distal_5                0                0
distal_6                0                0
distal_7                0                0
distal_8                0                0
distal_9                0                0
medial_0                1                0
medial_1                1                0
medial_2                1                0
medial_3                1                0
medial_4                1                0
medial_5                1                0
medial_6                1                0
medial_7                1                0
medial_8                1                0
medial_9                1                0
piece1_0                0                0
piece1_1                0                0
piece1_2                0                0
piece1_3                0                0
piece1_4                0                0
piece1_5                0                0
piece1_6                0                0
piece1_7                0                0
piece1_8                0                0
piece1_9                0                0
piece2_0                0                0
piece2_1                0                0
piece2_2                0                0
piece2_3                0                0
piece2_4                0                0
piece2_5                0                0
piece2_6                0                0
piece2_7                0                0
piece2_8                0                0
piece2_9                0                0
piece3_0                0                0
piece3_1                0                0
piece3_2                0                0
piece3_3                0                0
piece3_4                0                0
piece3_5                0                0
piece3_6                0                0
piece3_7                0                0
piece3_8                0                0
piece3_9                0                0
proximal_0              0                1
proximal_1              0                1
proximal_2              0                1
proximal_3              0                1
proximal_4              0                1
proximal_5              0                1
proximal_6              0                1
proximal_7              0                1
proximal_8              0                1
proximal_9              0                1
attr(,"assign")
[1] 0 1 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$mmVector
[1] "contr.treatment"

In [21]:
%%R -i minCounts

#keep <- filterByExpr(dds, min.count = minCounts, min.total.count = minCounts*5)
#dds <- dds[keep,,keep.lib.sizes=FALSE]
dim(dds)
[1] 25569    60
In [22]:
%%R
dds <- calcNormFactors(dds)
dds <- estimateGLMRobustDisp(dds, mm)
fit <- glmFit(dds, mm)

MAP contrasts!¶

In [23]:
levelsToMap
Out[23]:
['piece1', 'piece2', 'piece3', 'distal', 'medial', 'proximal']
In [24]:
%%R -o piece2_vs_piece1 -o piece3_vs_piece1  -o distal_vs_piece1  -o medial_vs_piece1  -o proximal_vs_piece1 -o proximal_vs_distal -o medial_vs_distal -o proximal_vs_medial
#-o Mid_vs_early -o Late_vs_mid
topn=25000

piece2_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=2), adjust.method = "fdr",  n = topn, sort.by = "logFC"))
piece3_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=3), adjust.method = "fdr", n = topn, sort.by = "logFC"))
distal_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=4), adjust.method = "fdr", n = topn, sort.by = "logFC"))
medial_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=5), adjust.method = "fdr",  n = topn, sort.by = "logFC"))
proximal_vs_piece1 = as.data.frame(topTags(glmLRT(fit, coef=6), adjust.method = "fdr",  n = topn, sort.by = "logFC"))
proximal_vs_distal = as.data.frame(topTags(glmLRT(fit, contrast=c(0,0,0,-1,0,1)), adjust.method = "fdr",  n = topn, sort.by = "logFC"))
medial_vs_distal = as.data.frame(topTags(glmLRT(fit, contrast=c(0,0,0,-1,1,0)), adjust.method = "fdr",  n = topn, sort.by = "logFC"))
proximal_vs_medial = as.data.frame(topTags(glmLRT(fit, contrast=c(0,0,0,0,-1,1)), adjust.method = "fdr",  n = topn, sort.by = "logFC"))

Volcanos.¶

In [25]:
DEGsDict = {}
DEGsDict = {"piece2_vs_piece1":piece2_vs_piece1,
            "piece3_vs_piece1":piece3_vs_piece1,
            "distal_vs_piece1":distal_vs_piece1,
            "medial_vs_piece1":medial_vs_piece1,
            "proximal_vs_piece1":proximal_vs_piece1,
            "proximal_vs_distal":proximal_vs_distal,
            "medial_vs_distal":medial_vs_distal,
            "proximal_vs_medial":proximal_vs_medial,
           }

LOGCFTHRESHOLD = 1.5
QTHRESH = 0.05
for k in list(DEGsDict.keys()):
    DEGsDict[k]["-logPVal"] = -np.log10(DEGsDict[k].PValue)
    DEGsDict[k]["significant"] = "notSignificant"
    DEGsDict[k].loc[(DEGsDict[k]["logFC"] < -LOGCFTHRESHOLD) &  (DEGsDict[k]["FDR"] < QTHRESH),"significant"] = "Down"
    DEGsDict[k].loc[(DEGsDict[k]["logFC"] > LOGCFTHRESHOLD) &  (DEGsDict[k]["FDR"] < QTHRESH),"significant"] = "Up"
    color_discrete_map = {'notSignificant': 'white', 'Up': 'red', 'Down': 'blue'}
    fig = px.scatter(DEGsDict[k], x="logFC", y="-logPVal",
                     #hover_data=DEGsDict[k][["genes","knownMarker"]], 
                     hover_data=DEGsDict[k][["genes"]], 
                     width=600, height=600, 
                     color = "significant", color_discrete_map=color_discrete_map,
                     #symbol="knownMarker",
                     title=k,
                    template="simple_white")
    
    
    fig.update_traces(mode='markers', marker_line_width=1, marker_size=10)
    
    fig.add_hline(y=DEGsDict[k].loc[DEGsDict[k]["FDR"] - 0.01 < 0,"-logPVal" ].min(), 
                  line_color="black",line_width=2, line_dash="dash")

    fig.add_vline(x=LOGCFTHRESHOLD, 
                  line_color="black",line_width=2, line_dash="dash")
        
    fig.add_vline(x=-LOGCFTHRESHOLD, 
                  line_color="black",line_width=2, line_dash="dash")
    fig.show()
In [26]:
for k in DEGsDict.keys():
    DEGsDict[k][DEGsDict[k].significant != "notSignificant" ].sort_values("logFC", ascending = False).to_csv(outdir+FinaLeaf+"/"+output_basename+"."+ k + ".tsv", sep = "\t")
    DEGsDict[k][DEGsDict[k].significant != "notSignificant" ].sort_values("logFC", ascending = False).to_excel(outdir+FinaLeaf+"/"+output_basename+"."+ k + ".xlsx")

Proximal and medial clinics¶

In [27]:
proximal_vs_distalDF = DEGsDict["proximal_vs_distal"]
proximal_vs_distalDF["Contrast"] = "proximal_vs_distal"

medial_vs_distalDF = DEGsDict["medial_vs_distal"]
medial_vs_distalDF["Contrast"] = "medial_vs_distal"

OverallContrsts = pd.concat([proximal_vs_distalDF,medial_vs_distalDF], ignore_index=True)
#OverallContrsts = OverallContrsts[OverallContrsts["significant"] != "notSignificant"]
gp = OverallContrsts.groupby(["genes","significant"], as_index=False).size()
UpAtl1 = gp[(gp["significant"] == "Up") & (gp["size"] >=1)].genes.tolist()
DownAtl1 = gp[(gp["significant"] == "Down") & (gp["size"] >=1)].genes.tolist()
OverallContrsts = OverallContrsts[OverallContrsts.genes.isin(DownAtl1+UpAtl1)]
# OverallContrsts = pd.concat([OverallContrsts,pd.DataFrame({"genes":OverallContrsts.genes.unique().tolist(),
#                                                            "logFC":0,"significant":"notSignificant",
#                                                            "Contrast":"distal"})])
OverallContrsts["significant"] = pd.Categorical(OverallContrsts.significant, categories = ['notSignificant','Up', 'Down'],ordered=True)
OverallContrsts["Contrast"] = pd.Categorical(OverallContrsts.Contrast, categories = ['proximal_vs_distal','medial_vs_distal','distal'],ordered=True)
OverallContrsts["ABSlogFC"] = abs(OverallContrsts["logFC"])

#OverallContrsts = OverallContrsts[OverallContrsts["significant"] != "notSignificant"]
In [28]:
#sns.catplot(data=OverallContrsts, x="Contrast", y="logFC", hue="logFC")
fig, ax = plt.subplots(1,1, figsize=(20,15))
sns.axes_style("ticks")
sns.scatterplot(data=OverallContrsts.sort_values(['Contrast','significant']), palette=color_discrete_map,
				x="Contrast", y="logFC",  size="-logPVal", hue="significant", ax=ax, linewidth=.2,sizes=(100, 700),
				edgecolor="black")
sns.lineplot(data=OverallContrsts, x="Contrast", y="logFC",units="genes", estimator=None, alpha=.05, color="grey", ax = ax, zorder=-1)
sns.despine()
ax.legend(bbox_to_anchor=(1.22, 0.8), prop={'size': 30})
# Plot labels
nGenes = 10
# Get pixels size
bbox = ax.get_window_extent().transformed(fig.dpi_scale_trans.inverted())
widthPXL, heightPXL = bbox.width, bbox.height
heightPXL *= fig.dpi


width = 1 
height = 1
QuadrantFractionDictH = dict(zip(["Down","notSignificant","Up"],
[1,0,-1]))



offset = lambda p: transforms.ScaledTranslation(p/72.,0, plt.gcf().dpi_scale_trans)
trans = plt.gca().transData

ContrastTuple = (OverallContrsts.loc[OverallContrsts['significant'] != "notSignificant",["Contrast",'significant']].drop_duplicates())
ContrastTuple = list(zip(ContrastTuple.Contrast, ContrastTuple.significant))
for i in ContrastTuple:
	PointsText = OverallContrsts[(OverallContrsts["Contrast"] == i[0]) & (OverallContrsts["significant"] == i[1])].sort_values("ABSlogFC", ascending=False).head(nGenes)
	PointsText["TextRank"] = (-PointsText["ABSlogFC"]).rank()-1
	PointsText["yText"] =  [QuadrantFractionDictH[i] for i in PointsText["significant"]]
	PointsText["yText"] = PointsText["TextRank"]*PointsText["yText"]
	PointsText = PointsText.sort_values("TextRank", ascending=True)
	PointsText["FCdelta"] = PointsText["logFC"] - PointsText["logFC"].values[0]
	PointsText["yText2"] = PointsText["yText"] - PointsText["FCdelta"]
	PointsText["xTextSwitch"] = np.where(PointsText["Contrast"] == "medial_vs_distal",-1,1)
	PointsText["xAnchorSwitch"] = np.where(PointsText["Contrast"] == "medial_vs_distal","right","left")
	PointsText["xAnchorSwitchArrow"] = np.where(PointsText["Contrast"] == "medial_vs_distal",1,0)
#point["xAnchorSwitch"]
	for i, point in PointsText.iterrows():
		ax.annotate(str(point['genes']), va="center", ha=point["xAnchorSwitch"],
					xy=(point["Contrast"], point['logFC']), 
					xytext=(200*point["xTextSwitch"], (point["yText2"]*heightPXL/26)),
					textcoords='offset points',fontsize=25,
			arrowprops=dict(facecolor='black', arrowstyle='-',relpos=(point["xAnchorSwitchArrow"],0.5)))




plt.title("DEGs: {} Qval:{} MinlogFC:{}".format(output_basename.split("_")[1], QTHRESH,LOGCFTHRESHOLD), fontsize=30)
plt.tick_params(axis='both',  labelsize=30)
plt.xlabel('Contrast',  fontsize=40)
plt.ylabel('logFC',  fontsize=40)

fig.savefig("./figures/Volcanos-{}_Qval-{}_MinlogFC-{}.png".format(output_basename.split("_")[1], QTHRESH,LOGCFTHRESHOLD), bbox_inches='tight',dpi = 400)
In [29]:
output_basename.split("_")[1]
Out[29]:
'Cycling'

number of DEGs¶

In [30]:
proximal_vs_distalDF = DEGsDict["proximal_vs_distal"]
proximal_vs_distalDF["Contrast"] = "proximal_vs_distal"

medial_vs_distalDF = DEGsDict["medial_vs_distal"]
medial_vs_distalDF["Contrast"] = "medial_vs_distal"

OverallContrsts = pd.concat([proximal_vs_distalDF,medial_vs_distalDF], ignore_index=True)
#OverallContrsts = OverallContrsts[OverallContrsts["significant"] != "notSignificant"]
# OverallContrsts = pd.concat([OverallContrsts,pd.DataFrame({"genes":OverallContrsts.genes.unique().tolist(),
#                                                            "logFC":0,"significant":"notSignificant",
#                                                            "Contrast":"distal"})])
OverallContrsts["significant"] = pd.Categorical(OverallContrsts.significant, categories = ['notSignificant','Up', 'Down'],ordered=True)
OverallContrsts["Contrast"] = pd.Categorical(OverallContrsts.Contrast, categories = ['proximal_vs_distal','medial_vs_distal','distal'],ordered=True)
OverallContrsts["ABSlogFC"] = abs(OverallContrsts["logFC"])






OverallContrstsSS = OverallContrsts[(OverallContrsts["significant"] != "notSignificant") & (OverallContrsts["Contrast"] != "distal")]
OverallContrstsSS["Contrast"] = OverallContrstsSS["Contrast"].astype("string")
OverallContrstsSS["significant"] = OverallContrstsSS["significant"].astype("string")
OverallContrstsSS = OverallContrstsSS.groupby(["Contrast","significant"], as_index=False).size()

OverallContrstsSS["Contrast"] = pd.Categorical(OverallContrstsSS.Contrast, categories = ['proximal_vs_distal','medial_vs_distal'],ordered=True)
In [31]:
import matplotlib.patches as mpatches

sns.axes_style("ticks",{'axes.grid' : False})
fig, (ax1, ax2) = plt.subplots(ncols=2, figsize=(20, 10), sharey=True, gridspec_kw={'wspace': 0})
sns.barplot(data=OverallContrstsSS[OverallContrstsSS['significant'] == 'Up'], x='size', y='Contrast',color="red",
            ci=False, orient='horizontal', dodge=True, ax=ax2)
ax2.yaxis.set_label_position('right')
ax2.tick_params(axis='y', labelright=True, right=True)
ax2.set_title('  '+'Upregulated Genes', fontsize=40)
ax2.set_xlabel('number of genes')
ax2.xaxis.get_label().set_fontsize(40)
ax2.yaxis.get_label().set_fontsize(40)

ax2.grid(False)
ax2.tick_params(labelright=False, right=False)
ax2.set_xlim(0,OverallContrstsSS["size"].max())
ax2.tick_params(labelsize=30)



sns.barplot(data=OverallContrstsSS[OverallContrstsSS['significant'] == 'Down'], x='size', y='Contrast',color="blue",
            ci=False, orient='horizontal', dodge=True, ax=ax1)
ax1.invert_xaxis()
ax1.set_title('  '+'Downregulated Genes', fontsize=40)
ax1.set(xlabel='number of genes')
ax1.xaxis.get_label().set_fontsize(40)
ax1.yaxis.get_label().set_fontsize(40)
ax1.tick_params(labelsize=30)

ax1.grid(False)
ax1.set_xlim(OverallContrstsSS["size"].max(),0)

sns.despine()
plt.legend(bbox_to_anchor= (1.6,.5), frameon=False, title="DEGs: {}".format(output_basename.split("_")[1]), fontsize=30,
              title_fontsize=30)

fig.savefig("./figures/NumberOfDEGs-{}_Qval-{}_MinlogFC-{}.png".format(output_basename.split("_")[1], QTHRESH,LOGCFTHRESHOLD), bbox_inches='tight')
No handles with labels found to put in legend.